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Abstract: A subalgebraic approximation algorithm is proposed to estimate from a set of time series the 
parameters of the observer representation of a discrete-time polynomial system without inputs which can 
generate an approximation of the observed time series. A major step of the algorithm is to construct a set 
of generators for the polynomial function from the past outputs to the future outputs. For this singular 
value decompositions and polynomial factorizations are used. An example is provided. 
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1. INTRODUCTION 

The system identification of polynomial systems is motivated 
by the need for models of biochemical reaction systems in the 
life sciences. Also in the area of control engineering and of 
economics there is a need to determine parameter values of such 
control systems from data. As far as the authors have been able 
to determine there is no satisfactory algorithm for the general 
problem of determining the parameter values of these systems. 

The problem of the paper is to determine a system in the class 
of discrete-time polynomial systems without inputs in the form 
of an observer realization such that it produces for each time 
series of outputs a predicted time series which is a reasonable 
approximation of the supplied output time series. 

The relevant literature on polynomial systems and their system 
identification is briefly summarized. At the time this paper is 
written there are available results on the realization theory of 
polynomial and of rational systems see Sontag (1979), Bar- 
tosiewicz (1988), Nemcova and van Schuppen (2009), 
Nemcova and van Schuppen (2010). The problem of structural 
identifiability of polynomial and of rational systems was solved 
by J. Nemcova in Nemcova (2010). Earlier papers of the au¬ 
thors include Nemcova and van Schuppen (2009), Nemcova 
et al. (2012). Various aspects of system identification of rational 
systems are also discussed in 
Boulier and Lemaire (2009), Gevers et al. (2013), 

Bazanella et al. (2014). 

The algorithm proposed in this paper determines a polynomial 
system in the form of an observer, thus driven by the available 
output. It will be proven using system theory that an observer 
polynomial system is equivalent with the conditions: (1) the 
state of the observer at any time is a polynomial function 
of the past outputs; (2) the future outputs are a polynomial 
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function only of the current state; and (3) the next state is a 
polynomial function only of the current state and the current 
output. The main step of the algorithm is to construct an 
approximate generator set for the polynomial equation from 
the past outputs to the future outputs. The subalgebra generated 
by this generator set is then an approximate subalgebra of the 
algebra of the function from past outputs to future outputs. 

2. PROBLEM FORMULATION 

System identification of a polynomial and rational system is 
motivated by the occurence of these systems in engineering 
(satellite orientation problems), biochemical reaction networks 
(mass action kinetics), and economics (products of prices and 
quantities). In this paper the authors focus on polynomial sys¬ 
tems. For the extension to rational systems and to systems with 
inputs there is insufficient space in this short paper. 

System identification is a research area that addresses the prob¬ 
lem of how to go from observational data to a system with 
its parameter values. The following procedure is often used: 
(1) Modeling. Model the phenomenon as a control system as 
understood in system theory; (2) Data collection. Collect data 
in the form of a time series of the phenomenon to be identified; 
(3) Structural identifiability. Determine whether the selected 
parametrization of the subclass of systems is structurally iden¬ 
tifiable and, if not, modify the system subclass such that the 
system subclass is structurally identifiable; (4) Approximation. 
Determine an algorithm to compute the parameters of a system 
in the considered subclass from the observation data which 
is a reasonable approximation according to an approximation 
criterion; (5) Complexity estimation. Determine that subclass of 
systems which achieves a reasonable value for the approxima¬ 
tion criterion and which minimizes the complexity. This paper 
addresses Step (4). 

For the approximation problem of system identification there 
are two main approaches: (1) Minimization of an approxima- 





tion criterion over the considered subclass of systems. This 
approach, though often used, suffers from the problem that 
the criterion is a nonconvex function of the parameters which 
makes the minimization approach practically almost impossi¬ 
ble. (2) An algebraic method based on realization theory of 
system theory. The method is explained in the next section. 

For the remainder of this paper the reader is expected to have 
read the notation and terminology of the appendices. 

Definition 2.1. A time series is a collection of real vectors 
denoted by its dimension d y £ Z+, its length in time steps 
fi £ Z+, and {y(t) £ R d », t £ {1,2, C Z}. Often 

there is available a finite set of time series. 

Definition 2.2. The system class considered is that of a discrete¬ 
time polynomial system without inputs in observer representa¬ 
tion, which produces y(t\t — 1), a prediction of the next output 
based on past outputs. 

x(t + 1) = f 0 {x(t),y{t )), z(0) = x 0 , 
y(t\t - 1) = h 0 (x(t)), 

T = {0,1,2,..., fi} c N, X = R n , Y = R d y, 
f 0 :X xY -+X, h 0 :X -+Y, 

where f Q and h a are both polynomial functions. 

Problem 2.3. Consider a finite set of time series and the sub¬ 
class of control systems of Def. 2.2. The problem is to deter¬ 
mine an observer system in the subclass considered, specified 
by its parameter values, such that the estimated observer sys¬ 
tem, when supplied with the output of the time series, produces 
a one-step prediction of the output which prediction time series 
is close to the observed time series. 

3. THE APPROACH 

The subalgebraic approach to system idenfication of polyno¬ 
mial systems is based on the subspace identification algorithm 
of Gaussian systems. The approach was initially proposed by 
H. Akaike based on contacts with R.E. Kalman, has been de¬ 
veloped for Gaussian systems, and is known as the subspace 
identification algorithm, see for references van Overschee and 
De Moor (1996). The subalgebraic approach to polynomial 
systems is based on realization theory of nonlinear systems in 
particular of bilinear, of polynomial and of rational systems. 
See for references on realization theory those mentioned in the 
previous section and Fliess (1990), Isidori (1973), Vidal (2008). 
No approximation criterion is used in the paper. The principle 
of the subalgebraic method is explained with the following 
theorem. 

Theorem 3.1. Consider a discrete-time polynomial system 
without inputs in its observer realization with as output the one- 
step prediction, 

x a {t + 1) = f 0 {x 0 (t),y(t)), x o (0) = x ofi , (1) 

y(t\t - 1) = h 0 (x 0 (t)), f 0 , h a , polynomial maps. (2) 

The observer system representation (1,2) may be transformed 
to the following set of polynomial functions assuming that 
the observer is a true observer, hence the predictions equal 

y(t) = y{t\t — 1), and in terms of ( y + {t), y~(t )) as defined in 
equation (8), 


x 0 (t) = gio(x ofi ,y~(t)), (3) 

y + (t) = h i0 (x 0 (t)) = h io (g io (x ofil y~(t))), (4) 

x 0 (t + 1) = f 0 (x 0 (t),y(t)), f 0 , gto, h io , polyn. (5) 

Note that by the equations (3,4) the current state is a polynomial 
of the components of the past outputs and the future outputs are 
polynomial functions of the components of the current state; 
and by equation (5) the next state x Q (t + 1) is a polynomial 
map in the tuple (x 0 (t),y(t)) of the current state and the current 
output. If one considers the initial state as a constant then the 
equations (3,4) define for any time a polynomial function from 
the past outputs to the future outputs. 

The algorithm for the subalgebraic approximation is based on 
the above theorem and consists of the steps: 

(1) Compute a state vector x(t) in terms of a polynomial 
function of past outputs such that the future outputs are 
a polynomial function of it. In terms of formulas, 

V + {t ) ~ fio(y~(t)) = h io (g io (y~(t))) = hi 0 (x(t)), (6) 
x(t) = gio(y~(t)), fio, 9io, hio polynomial functions. 

(2) Compute the polynomial system dynamics 

x(t + 1) « f 0 (x(t), y(t)). (7) 

The main task of the algorithm is to compute a set of generators 
of the polynomial map f lo from past outputs to future outputs. 
The complexity of the computations is limited by several steps 
of the algorithm. 


4. THE ALGORITHM 


Definition 4.1. The subalgebraic approximation algorithm for 
system identification of discrete-time polynomial systems. 
Data: A time series of outputs with the notations: the dimension 
of the output d y £ Z+, the length of the time series t\ £ Z+, 
the number of time series s £ Z+, and finally the time series 
matrix Y ts £ R* 1 x dy x s . The parameters of the algorithm are 
ri, r 2 , r 3 , r 4 £ (0,1) C R and (f+ in , (f+ ax , t ~ ax ) £ 

Z5_ defined below, and various maximal power vectors. 


(1) Construct the vectors of the past and of the future time 
series. Take a time t £ T = {1,2,..., G} less or equal 
to ti/2. Denote the length of the tuple of the future and of 
the past output time series respectively by ( t + , t~) £ Z + 
and set their extrema such that (t + i+ ax ), (i — f“ ax ) G T. 

Iterate from Step 2 to Step 6 in a Levinson-like manner 
by increasing the horizon lengths (t + ,t~) from the values 

(*min> *min) t0 the ValueS (4ax, *max) G Z +- 

Construct the symbolic vectors of the future and the 
past series, and their values for the each of the time series, 


d y r t tly, d y — t dy £ *3. , 

(y + (t),y~(t)) £ (R d » + xrV) 

(fv(t + t + - i)\ fy(t- l) \\ 

/ / ../X , J.+ <1\ 1 / ,„/+ O'! I \ 


( 8 ) 


y(t + f + — 2) 




y(t - 2 ) 


J V y(t — t )/ / 


(: y + {t,k),y~{t,k )) £ (RV x RV), k £ Z a . 

(2) Define the power matrices of the future and past output 
time series. 



k 


max,y 


d v + = d y +, K v + = Id v+ G N d « +Xd y + , 
^max,t/ eN d *,k* = max k max,y (i), (9) 

*€ Z dy 

( k T k T ) T G N d y~ 

V A 'max,y> • * * ’ "'max,!/ / e 5 


d y 

dv~ < ((^max,y(^) + 1))* < (ky + l) dy ~ , 

i=l 

K v -e N d «- x V(fc maXil ,_), 

(Z/„+ , AT„+ ) = {Id v+ 1 Id v+ )> 

(L„-, A'„-) = (7d <) _ , K v ~). 

(3) Iteration of blocks of the power matrix K v ~. If the row 
dimension of the bounded power matrix K v - of the past 
outputs is very high, say larger than 500, then execute 
this step. Partition the full power matrix K v - into a 
finite number of row blocks 7C_ ,..., Kff 1 ’ starting 
with lowest power blocks. Denote the corresponding row 
dimensions by d (1 J,..., d^ mi \ 

Iterate from Step 4 to Step 5. Start with the first block 
(cf 1 }. L ^}, K (l j). After each cycle add to the current 

generator set indexed by (iS 9 !, K ^ 9 }) the next block of 
the power matrix. 

(4) Construct the monomial vectors of the future and the past 
outputs. Construct next for the parameters 

(d v +, K v +, d^™\ K™) set in Step 2 or Step 3, the asso¬ 
ciated monomial vectors according to Def. B.2, 

v+(k) =v(y + (t,k),d y+ ,K v +) = y + (k ) G R d -+(10) 

V + (t) = (u + (l) v + (2) ... v + (s )) € R d "+ xs , (11) 

v-(k) = v(y-(t,k),d y -,K ( u m >) GRV , (12) 

d y~ 

' l ’i”( fc ) : =ri yj(t,k) K v- M , see (B.l), (13) 

j=i 

V~(t) = (v~{ 1) v~{2) ... v~(s)) el-’" 5 . (14) 

Note that the dimension and hence the complexity of v~ 
is exponential in terms of d v - = t~d y . 

(5) Reduce the generator set by (1) linear dependence. Com¬ 
pute, according to Algorithm E.4, the approximate mono¬ 
mial equation of the future and the past time series. 

(«i, D ni , C v +, Li,X, H*,tablei) (15) 

= SVDtnmc(dy + , dy~, dv y +, dv y -, s, 

V+(t),V-(t),n ); H* 

Tr ( TTX ) 

V + (t)*H*(t)V-(t)=C v+ L 1 (y-(t,*)fv- ,(16) 

C v + G R d "+ Xni , Lx € . (17) 

(6) LK-Reduction. Reduce the generator set further for the 
matrix tuple ( L \, K y _) by deleting those columns of the 
matrix L i and the corresponding rows of the matrix K \ 
whose Zi-norm of the column is less than G (0,1) times 
the 1-norm of L\. Thus delete column j of L \ if, 

n n 

^2\Lx{i,j)\ <r*||Li|| Ll = r * max^ \Lx(i,j)\, 

i=i 3&n i=l 

result. (18) 


(7) Reduce the generator set by (2) elimination of prod¬ 
ucts of generators. Compute a new generator set with 
possibly less generators according to Step 2 of the al¬ 
gorithm of Def. E.2. Starting from equation (18). with 

(C v +, (l[ 9 },k1 9 J)). The result is, 

y + (t) =v + {t) w C v +L { v 9 }(y-(t)) K v- 

~hf(x(t )) = L h +x(t) Kh to, (B.l) (19) 

x{t) =g io (v~(t)) = L gio y~ (t) Kgi ° G R™. (20) 

(8) Compute the output equation. 


y{f)=P v{t) y + {f) « Py(,)hi 0 (x(t )) = h 0 (x{t)), 

P y (t) € M. dvXd v + , a projection. 

The next steps aim at the computation of the system 
dynamics, see equation (7). 

(9) Compute the value of the next state. First compute the 
past time series at the next time index (t + 1). Secondly, 
compute the value of X(t + 1) for each time series. 


y + 


/ y(t,j ) 

I y(t-i,j) 


\ 

e R d y ~, 


( 21 ) 


\y(t-t~ + i ,j)) 

v~{t + l,j) =v(y-(t + hj),d y -,Kg io ) € K d o- , (22) 

V~(t + 1) = ( v~(t+ 1,1) v(t + 1,2) ... v(t + 1, s) ) (23) 
G R d v- Xs ; d Vx = n e Z+, 


i4(i + i) = x(i + i) = L Jio r(( + i)er Xs . 

(10) Compute the monomial vector of the current state and the 
current output. 


d (x, y ) = n + d y , d V(x y) € Z+, (24) 

n dy 

dv( x y) = JJ[fcmax,a:(») + 1] [fcmax y2 (j) + !]> ( 2 5) 

1=1 j= 1 

choose (max,! G N , fc m ax,y2 G N 

^ma x,(x,y) = ( ^max,x^max,y2 ) ^ ^ ( X ’V) , 

Rv (x _ v) G Xd(x ’ v) (k max , (x ,y)), 

v (x(t),y(t))(h) = v ((y(t’kj) >( n + dy),K V(xy) ) (26) 

V *,v = { v ( X (t),y( t ))( 1) ••• V (x(.t),y(t)){s) ) ■ ( 27 ) 

Note that the complexity of the expression of d Vx is 
exponential in n + d y . 

(11) Reduce the generator set by (1) linear dependence. Com¬ 
pute the approximate polynomial function of the next fu¬ 
ture state depending on the vector of the current state and 
of the current output. Compute according to Def. E.4, 

(n 2 , D n2 ,C 2 ,L 2 , X 2 , P 2 , table 2 ) (28) 

= SVDtrunc(n, d ( x . y ), n, dv ^ v ), s, V x (t + 1), 

V x , y ,r 2 ); L 3 = H2 GR BXd ’c..). 

X(t+1) = V x(t+1) w 7'3^(x(t),2/(t))| (29) 

x(t+l) = L 3 (x{t),y(t)) K3 . (30) 

(12) Approximate the monomial map. Reduce the matrix tuple 
(L 3 , K 3 ) to ( Lf io , Kf io ) as in Step 6 for the matrix pair 

{L 2 ,K% 



(13) Compute the polynomial system. Finally compute the 
discrete-time polynomial system without input in the ob¬ 
server form, by writing the linear map of nomomials as a 
vectorial polynomial function, 

f 0 (x,y) = Lf io (x,y) Kf *o, see(B.l), (31) 

x a {t + 1) = fo(x 0 (t), y{t)), £ o (0) = x 0 ,o, (32) 

y Q {t\t - 1) = h 0 {x 0 (t)), from Step (8), (33) 

X ofi = X(t) G R nxs . (34) 

Output (n, f 0 , h 0 , {X ofi (k), y + {k), y~{k), k G Z s j). 

The user is advised to take t + and t~ equal to about 4 times the 
dimension of the expected state space. 

5. EXAMPLES 


Example 5.1. For a set of times series generated by a simple 
polynomial system with a one-dimensional output and a two- 
dimensional state vector, a computer program for Algorithm 4.1 
has computed the following observer polynomial system. 


x(t + 1 ) = f 0 {x(t),y{t)) = L fo {x 0 {t),y{t)) K t°, 


y{t\t - 1) = h 0 (x(t)) = Cx 0 {t), C = ( -0.0225 0.0336), 


( 0.009 .089 0.023 .571 -.004 -0.020\ 
( —0.015 .309 —0.014 .074 .008 0.212 ) 


(x 0 , y) Kf ™ = (x 0 ,ix 0t2 y, x 0 px 0 , 2 , x 0 py, x Q p,x a py, x 0 p) T . 
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Appendix A. NOTATION 


The set of the integers is denoted by Z and the positive integers 
by Z + = {1,2,...}. The set of the natural numbers is denoted 
by N = {0,1,2,...} and by N™ its ?r-tuples. For any n G Z + , 
denote Z„ = {1,2,..., rz}. The set of the real numbers is 
denoted by R, the positive real numbers by R + = [0, oo), and 
the strictly-positive real numbers by R s+ = (0, oo). The vector 
space of n-tuples of the real numbers is denoted by R n , for 
n G Z + . The set of matrices with entries in the real numbers of 
size n x to, for n, m G Z + , is denoted by R raxm . A diagonal 
matrix of the set of square real matrices R" xn for n G Z + is 
a matrix such that its off-diagonal elements all zero and the set 
of such matrices is denoted by R)} x ™. The subset of positive 
diagonal matrices is defined by the condition that D,p > 0 for 
all i G Z„ and it is denoted by K-^iag +• Similarly, K-^iag s + - 

Definition A.l. The truncation operation of a positive diagonal 
matrix based on the i \-norm of the diagonal. 

Data. (n,D,r) G (Z + x Rd;*g i+ x (0) 1))> D ^ 0, where r 
is an approximation threshold. Assume that Dip > I) 2:2 > 
... > > 0 . 

(1) Compute the /1 -norm of the diagonal elements of the 
diagonal matrix D, ||diag(D)|| il = E"=i D l/t . 

(2) Compute 

n r = argmin ieZn {X;- = i A,*/ ||diag(D)|| Zl > r}. 

(3) Construct the approximant positive diagonal matrix 


On r ,i,i — Dip,Vi G {1,2,..., n r f, D n 


diag,+ ’ 


D r 


D nr 0\ 

0 0 J 


G R 


nxn 
diag, +' 


(A.l) 


(4) Output (n r , D r , tablei) G (N x R^}} + x R rax2 ), where 
tablei = {(j,Ei=i A,i/||diag(£))||j 1 ), j G Z„}. 


Appendix B. MONOMIALS AND MONOMIAL VECTORS 

Definition B.l. Consider a set of n G Z + commutative vari¬ 
ables denoted by x = (x\,... ,x n ). A monomial is a term of a 
polynomial defined by the formulas. 



= n< W =*f ( 1 ) * 2 ( 2 ) ---*n (n) . 




G m onN = {x k \V k £ N"}, A(G mon )M = A R „ [®], 
p{x) = J2 < k ) xk e R M> v k e N", c(fe) € R. 


fceN" 


Call a:* a monomial in the indeterminates x, a vector /c £ N" 
a. power vector , N" a set of power vectors, G mon M the set of 
all monomials of x, and p £ R [a:] a polynomial in monomial 
representation. 

There exists a bijective correspondence between the set of 
power vectors N" and the set G mon M of monomials in the 
indeterminates x, described by the map k x k . The partially- 
ordered set N" may be equiped with a monomial ordering. 
Below the specific monomial ordering called the lexicographic 
order relation on the index set of power vectors N" will 
be used, see (Cox et al., 1992, Def. 2.2.3). It is denoted 
by >iex- By the bijective correspondence between N" and 
GmonM the lexicographic ordering of N" is transformed into 
a lexicographic ordering on G rn0 n [a-'] which is again denoted 
by >iex and which will be called the lexicographic ordering of 
G m on[at]. 

Definition B.2. Define the power-bounded monomial vector of 
a set of n £ Z + commutative variables by the following 
formulas. 


h cz M" 

ttmax vz J-n , 

N n (k max ) = {k £ N"|0 < k{i) < k max (i), v i £ Z n }, 

n 

k* = |N n (fc max )| = n( fc max(i) + 1) € Z+, 

2=1 


choose d v G Z+, 0 < d v < k*, 
N d,,xn (fc max ) = {K £ N d ” xra |0 < K(i,j) < ft max (j)}, 
choose K v £ N d,,xn (fc max ), and define, 
v = v(x, n, K v ) = x Kv £ R d ”, (B.l) 

v(x,n,K v )i = xf v ^’^x^ v ^’ 2 ^.. .x^ v ^ ,n ^; h £ R d ” , 

p(x) = c(k)x k = h T v(x,n, K v ) = h T x Kv . 

keK v 


Call fc max the maximal power vector and fc m ax(*) the maximal 
power of xp, N™(fc max ) the bounded power vector set , which 
set inherits the lexicographic ordering of the elements of N n ; 
K v £ N d ” xra (fc max ) the power matrix of the monomial vec¬ 
tor v(x, n, K v ), where the row elements of K v are the pow¬ 
ers vectors of the vector v in decreasing lexicographic order; 
v(x,n, K v ) a (symbolic) monomial vector, which contains all 
monomials indexed by K v £ N d,,x ™(fc max ) in their lexico¬ 
graphic order from high to low order; the number of compo¬ 
nents of this monomial vector equals d v < k*\ and finally 
call the equation p(x) = h T v(x, n, K v ) = h T x Kv , the power- 
bounded monomial representation of the polynomial p. 

Example B.3. There follows an example of a monomial vector. 


x = (xi, x 2 ), n = 2, k max — (2,1) T , 


v(a,2,(2,l)) 


dy = k* = (2 + 1) x (1 + 1) = 6, 


/ x\x 2 \ 
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1\ 

xl 


2 

0 

X\X 2 

• K — 

1 
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Xl 
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Appendix C. POLYNOMIALS 


In this appendix and the next one several aspects of the commu¬ 
tative algebra of polynomials maps are described including al¬ 
gebraic geometry. See Becker and Weispfenning (1993a); Moor 
(2014), Cox et al. (1992), and Zariski and Samuel (1958). 

Let n £ Z + be a positive integer. The ring of polynomials in 
n variables with real coefficients is denoted by Rjaq,..., x n \. 
The simplified notation of R[x] will be used if it is understood 
that the variable x has n components. Examples are 2a: 2 + 3a: + 
4 £ R[x] and 2lx\x 2 + llx\x 2 + lx 2 £ R[a;i, x 2 \. 

Below polynomial functions of tuples of the real numbers 
X = R™ are needed. A polynomial function on R” is a map 
p : X —> R n for which there exists a set of polynomials 
qi,,q n £ R[aq, x 2 ,... , x n ] such that Pi = (p on R” for all 
i £ Z n . Denote by Ax the set of all polynomials on X = R™. 

Definition C.l. The monomial representation of a finite set 
of polynomials in the indeterminates x = x±... x,i r power 
bounded by the vector fc max £ N n , where the polynomials are 
the components of Lx K , is defined by the notation. 


G = 


Lx 


K 


i,card G xd„ 


Ml (l,k) 


xN d ” x4 (U)) 


Example C.2. Consider the set of polynomials, 

, K = 


L = 


0.1 0.2 
0.3 0.4 


3 1 
1 2 


K _ ( x\x 2 


(C.l) 


\Xix\ ) 

n(r ) = Tr K = ( O.latfats + 0.2a:ia:^ \ 1 

V 0.3xlx 2 + 0.4a:ia:2 J J ' 


Appendix D. A SET OF GENERATORS 

A subalgebra A\ of the algebra Ax is a subset A\ C Ax such 
that the algebraic operations of A\ are obtained from those of 
Ax by restriction and such that it is also an algebra itself in 
terms of those operations. 

Consider a subset G C Ax- The smallest subalgebra of A\ 
containing G exists, it is called the algebra generated by the 
set G, and it is denoted by Ax(G) C A\. A subalgebra 
Ai C Ax is called finitely generated if there exists a finite 
subset Gf C Ax such that A\ = Ax{Gf). 

Definition D.l. Call the finite set G, see (C.l), a generator 
set of the subalgebra A C R[x] and any row component of 
Lx k £ G a generator of A if A = A{G). Call it a nontrivial 
generator set if no column of the matrix L in the representation 
is entirely zero. Call G a minimal generator set of the algebra 
A = A(G) if it is nontrival and for any other generator set H it 
holds that cardc < card h- A set of generators of a finite set of 
polynomials H is a finite set of polynomials G such that 



A{H) = A(G), where, H = {p u ... ,P car( i ; , e K[x]}, 

G = {<?i, • ■ • iflcardc e R WI- 

The set of real numbers is also a ring. A finite subset 
{pi,... ,Pk} C R[x] is called algebraically dependent over 
R if there exists a nonzero polynomial / £ R[p] such that 
f(pi, ■ ■ ■ ,pk ) = 0. It is called transcedental otherwise. See 
(Zariski and Samuel, 1958,1, §17, p.28). An extension not de¬ 
scribed here because of lack of space is to define a transcedence 
basis for the algebraic structure used which then allows the use 
of the algorithms of Mtiller-Quade and Steinwandt (2000) for 
the computation of such a basis. 

Problem D.2. Consider a finite set of polynomials H C R[x], 
Construct a minimal set of generators G C t[x] of H. 

The above problem is equivalent to the problem of constructing 
a polynomial factorization of a polynomial map as defined next. 
Definition D.3. A polynomial factorization of a polynomial 
map y = f(u) with d y , d u £ Z+, / : R d “ —> R^, is defined 
to be a factorization of the form, 

y = /(«) = Hg( u )) = H x )> h e R[x], (D.i) 

x = g{u) acR dx , g£ R[u| (D.2) 

Gf = {h,...,fdy £%]}, Gg = {gi,...,gd x £ »[«]}, 
A(G f ) = A(G g ), (D.3) 

hence G g is a set of generators of A(Gf). Note that G g is 
a minimal set of generators if d x £ Z + is minimal over all 
factorizations. 

Appendix E. APPROXIMATION OF A GENERATOR SET 

Problem E.l. The problem of approximate polynomial factor¬ 
ization. Consider a polynomial function y = f(u) as defined 
in Section C. (The notation of this section differs from the 
main body of the paper.) Determine an approximate polynomial 
factorization of the form, 

V = f(u) « h(g(u)) = h(x), x = g(u), (E.l) 

d x £ Z + , g £ R[u], h £ R[x], 

The approximation criterion of the expression [ f(u ) — h(g(u))] 
is not specified. 

Definition E.2. The approximate polynomial factorization al¬ 
gorithm. Data, y = f(u), d y , d u £ Z+, / £ R[w], 

(1) Construct the approximation consisting of a linear-poly¬ 
nomial factorization by the algorithm of Def. E.4, 

V = fW) ~ C r g r (u) = C r x r , x r = g r (u), (E.2) 
d Xr £ Z + , x r £ R dxr , C r £ R d « xd *r. (E.3) 

(2) For the linear-polynomial factorization of Step 1 construct 
a new generator set of, possibly, lower cardinality than 
before, or, equivalently, a polynomial factorization, 

y = /(«) ~ C r g r (u ) = h(g(u)) = h(x), (E.4) 

x = g{u), (E.5) 

d x £ Z+, x £ R dx , h £ R[x], g £ R[u], 

by polynomial factorizations of the components of g r . For 
example, if for polynomial g r/n , there exists a factorization 
of the form g r ,m — . 7 r.i g r.j A 9r,s> where fjr.i - gr,j ? 9r,s 


have lower powers than those of g r , m , see (Becker and 
Weispfenning, 1993b, Sec. 5.1.). 

Comments follow. (1) A Grobner basis algorithm is not appro¬ 
priate for Problem E.l of the polynomial function y = /(it) 
because that would not be an approximation. In addition, the 
computational complexity is too high. (2) A Grobner basis 
algorithm may well be useful for Step 2 of Algorithm E.2. 
For the computations a simple procedure is used for Step 2 of 
Algorithm E.2. The procedure is not described here in detail 
because of lack of space. 

The transformation for a polynomial factorization is briefly 
described. Note that if 


V = f{u) ~ C r g r (u), x = g r (u), (E. 6 ) 

9r,m(fi) — x r,m — x r,i x r,j = 9r,i(fi)9r,jifi) 5 then, (E.7) 

X r — ( Xi . . . X m —i X m -l-l • ■ • X d x ) 1 

X — ( Xi . . . X m —i Xj.Xj X m -\-l ■ * * X d x ) > 

= Px ( XjXj Xi . . . Xd x ) = P x v{xfj = P x Xr 5 

y « C r x = C r P x v{x r ) = Cv(x r ) = h(x r ), (E.8) 
g(u) = Pg r g r (u), (E.9) 

y ~ C r x = C r g r (u ) = h(x r ), x r = g(u), (E.10) 


It is not claimed that the above procedure determines a minimal 
generator set which is not true in general. 

Definition E.3. The monomial equation of data matrices. Con¬ 
sider the polynomial equation y = f(u) and its monomial rep¬ 
resentation. Consider the case in which one is provided several 
tuples of values of input and output vectors, {{fjj , Ui) £ Y x 
U\ i = 1, 2,..., s}. A monomial equation of the data matrices 
for this set of tuples is then a linear map represented by the 
coefficient matrix H. 

Vy = HV U , H £ R^w , (E.l 1) 

V v = (v{y 1 ,d y ,K Vv ) ... v(y s ,d v ,K Vy )) eR^' 5 , 

V u = (v(ui,d u ,K Vu ) ... v(u s , d u , K Vu )) £ R d ”“ xs . 

Definition E.4. Linear approximation of a polynomial map. 
(Golub and Loan, 1983, Sec. 6.1). 

This algorithm is called SVDtrunction in Steps 5 and 11 of 
Algorithm 4.1. Data. (d y ,d u ,d Vy ,d Vu ,s,V y ,V u ,r) £ (Z+ x 
R d ”t/ Xs xR d ”“ xs x (0,1)). 

(1) Compute the singular value decomposition of the data 
matrix of the inputs. 


V u = vf SV 2 £R dv " 

a xs 

1 

(E.l 2) 

V\ £ R^” 1 2 * xdvv 

Vi £ R sxs , orthogonal matrices. 

s =(n) eRj 

”“ xs , ni gN, 


D = diag(di, d 2? • • 

d\ ^ ^ • • • 

., d ni ) £ 

> d ni > 0. 



(2) Compute according to the algorithm of Def. A. 1, of the di¬ 
agonal matrix D its truncation D n upto the approximation 
fraction. 



(E.13) 


(n, D n , tablei) = mdtrunc(ni, D), 

S '+ = ( ^ 1 o ) € K sxd -, D n G R^ ag " s+ , 

v+ = v. 2 T s+v \ g R sxd ““ , zr = v„v+ g R d ”« xd «« . 

(3) Compute the factorization according to, 

H* = [V y V 2 T ( 1 )][(/„ 0) 14] = CL, (E.14) 

L = ( I n 0) Vl G R nxd ”“ , X = LV U G R nxs , 

C = V v V.J (^Q nl ^j eR d ”« xn , (E.15) 

Vj, &H*V U = CLV u = CX. (E.16) 

(4) Output (n, D n , C, L , X, fP, table) with 
table = {(j,£i=i D M /||diag(D)|| il ), j GZ ni }. 



